Predicting compatibilized polymer blend toughness

Polymer blends can yield superior materials by merging the unique properties of their components. However, these mixtures often phase separate, leading to brittleness. While compatibilizers can toughen these blends, their vast design space makes optimization difficult. Here, we develop a model to predict the toughness of compatibilized glassy polymer mixtures. This theory reveals that compatibilizers increase blend toughness by creating molecular bridges that stitch the interface together. We validate this theory by directly comparing its predictions to extensive molecular dynamics simulations in which we vary polymer incompatibility, chain stiffness, compatibilizer areal density, and blockiness of copolymer compatibilizers. We then parameterize the model using self-consistent field theory and confirm its ability to make predictions for practical applications through comparison with simulations and experiments. These results suggest that the theory can optimize compatibilizer design for industrial glassy polymer blends in silico while providing microscopic insight, allowing for the development of next-generation mixtures.


SELF-CONSISTENT FIELD THEORY (SCFT) PARAMETERIZATION
For the SCFT parameterization in the model to hold, four assumptions must be met.One is system-dependent (assumption 1), and another has solid theoretical underpinnings (assumption 2).The other two assumptions require further justification.
We first check assumption 3, that the entanglement length (N e ) does not change as a function of distance from the interface (z) in Fig. S1A.This plot shows N e (z) for a selection of compatibilizer block lengths (L) for our simulations with an areal density of compatibilizers Σ = 0.045σ L − J 2 and Flory-Huggins parameter χN = 650.Here, σ LJ is the length scale of the Lennard-Jones minimum.Interestingly, it shows positive and negative fluctuations away from the bulk value of N e depending on the L and z.While this image clarifies that assumption 3 does not hold exactly, fluctuations away from N e are less than 30% in all cases that are not immediately adjacent to the interface and often are short ranged.These deviations are even less at lower Σ.Thus, we believe that N e (z) = N e is a reasonable assumption in many circumstances.
We next move onto assumption 4, which is entanglements' positions are independent of their position along the chains.To test this claim, we consider the normalized correlation of entanglements a distance s from a given entanglement on a compatibilizer, P e,C (s) is the probability of monomer s of a compatibilizer entangling, and the average is taken when an entanglement occurs at s ′ .This function takes a value of 1 when entanglements are completely correlated and 0 when entanglements are completely decorrelated.As shown for a representative compatibilizer in Fig. S1B, we see no correlations other than self-correlations (s = 0) in our data, indicating that assumption 4 holds in all cases.
To understand how much assumption 2 affects the model's predictions, we plot the true and predicted microscopic quantities in Fig. S2.First, we consider how these errors affect the AB crazing energy by considering the total length of HCH load-bearing strands per area, M HCH ⟨N ⟩ HCH , in Fig. S2A.The SCFT parameterization overpredicts these values by around 20%.As this term consists of two parameters M HCH and ⟨N ⟩ HCH , we can further ask which of these variables causes this overprediction.By plotting the true and predicted values of M HCH , Fig. S2B demonstrates that M HCH is the source of these issues.Thus, we can conclude that the SCFT parameterization overpredicts AB crazing energy because it predicts too many HCH load-bearing strands.Next, we move to the AC-crazing energy.We show that the SCFT model underpredicts the number of AC entanglements in Fig. S2C by about 10%.Thus, the parameterization underpredicts AC crazing energy.

DIBLOCK COMPATIBILIZERS
Diblock copolymers represent the simplest possible compatibilizer.As such, they are a valuable case study.First, we consider the scaling of toughness as a function of the areal density of compatibilizers (Σ) for long diblock compatibilizers at χN = 650 with a reference volume of the size of the monomer (= 1/ρ A = 1/ρ B , where ρ m is the number density of species m) in Fig. S3 to determine if it is consistent with Brown's theory which predicts toughness scales as Σ 2 (28).The solid guiding line demonstrates that our data aligns with this scaling.Next, we consider the SCFT parameterized theory.This model goes to zero for the smallest Σ but aligns closely with Brown's predictions and the data in the other cases.
We now turn to whether one needs to consider load-bearing entanglements rather than all entanglements.Here, we consider short diblock compatibilizers with a block length of L = 80 monomers and homopolymers of length N = 256 monomers.We fix N e = 10.3 monomers, set χ N = 650, and vary the areal density of these compatibilizers from Σ = 0.036/σ 2 to Σ = 0.18/σ 2 .As in Fig. 3  making them an excellent test of whether to keep the first and last primitive paths.When we use the original model with load-bearing entanglements, the true toughness values are well-predicted by the model that works with the significantly longer compatibilizers, as shown in Fig. S4.Because these compatibilizers are not so affected by removing their first and last primitive paths, this matching indicates that using load-bearing entanglements in our model is reasonable.On the other hand, when we apply the model with all entanglements, toughness is significantly overpredicted in all cases.Thus, using load-bearing entanglements allows our theory to work for considerably smaller compatibilizers.

TOUGHNESS OF POLYMERS OF THE SAME LENGTH
While Fig. 4C increases N such that N/N e is roughly constant, one may ask how toughness changes when N is kept constant.This question is not straightforward.In Fig. S5, we perform simulations in which we vary N e by varying chain stiffness but keep polymer length fixed at N = 256 monomers for two different block lengths L = 128 and L = 32 monomers.Figure S5A shows that if the block-copolymer blocks are long enough (L ≫ N e ), blend toughness increases with N e because this will lengthen the load-bearing strands.However, Fig. S5B illustrates if the blocks are short, changing N e can cause a decrease in the number of load-bearing strands.This effect occurs because the blocks are not long enough to entangle with the homopolymer matrix, causing blend toughness to decrease.

MONTE CARLO PROCEDURE FOR INITIAL POLYMER CONFIGURATIONS
For our initial homopolymer configurations, we fix bond length to r b , the minimum of the FENE potential.Starting with a single monomer, we grow polymers of length N by adding monomers sequentially such that they have a bond length r b and a polar bond angle (θ) chosen from the distribution where T = ϵ LJ /k B is our equilibration temperature, and l (θ) = r b 2 + 2 cos (θ).Here, ϵ LJ and σ LJ are the well-depth and length scale of the Lennard-Jones potential in our molecular dynamics simulations.The first term of Eq.S4 controls the polymer's stiffness, while the second term of Eq.S4 approximates the non-bonded interactions of the next-nearest neighbor, making our random walk non-reversible.Given a fixed r b , the distance between a monomer and its nextnearest neighbor is l (θ).We fix the stiffness of these polymers, κ, to match the end-to-end distances found in Auhl et al.We begin the initialization of our compatibilizers using the same method as the homopolymers, treating all monomers as the same type.We then place their center of mass on the interface (z = 0).To equilibrate them fully, we perform Monte Carlo moves on our old compatibilizer configurations to generate new configurations.The Monte Carlo moves we choose are (i) solid rotations, (ii) azimuthal rotations of bond angles (ϕ), and (iii) polar rotations of bond angles (θ).The solid rotations of the compatibilizer from the initial configuration are chosen from a Gaussian distribution with a standard deviation of 0.05 radians.Azimuthal rotations are chosen from a Gaussian distribution with a standard deviation of 0.35 radians from the initial angle.Angles for the polar rotations of bond angles are redrawn from Eq. S3.After each of these moves, the center of mass is shifted back to the interface.
Assuming a homogeneous field of homopolymer A for z < 0 and homopolymer B for z > 0, it can be shown that the energy of an A monomer a distance z from the interface is where U 0 = 2πϵ LJ ρσ 3 LJ is an overall energy scale, ϵ AB LJ is the Lennard-Jones well-depth of the A-B species, and ρ = ρ A = ρ B is the system density of species A and B. We perform this derivation by integrating the non-bonded Lennard-Jones interaction over all space except a sphere of excluded volume around the monomer with a σ LJ radius.A similar expression can be derived for a B monomer, U B (z), where the second term has the opposite sign.Thus, the non-bonded energy compatibilizer energy is where m i is the type of monomer i and z i is the distance of monomer i away from the interface.The Monte Carlo acceptance criteria is then We repeat these Monte Carlo moves 10 6 times for each compatibilizer, at which point U n.b. plateaus.We randomly place these compatibilizers in the xy plane at the interface.As a final step in generating our initial configurations, we perform 300 small solid rotations and translations to each homopolymer to minimize density fluctuations in the box as suggested in Auhl et al. (55).

FLORY-HUGGINS PARAMETERS OF MOLECULAR DYNAMICS SIMULATIONS
We run additional simulations of our polymer blends without compatibilizer to obtain our molecular dynamics simulations' Flory-Huggins parameters.We measure the interfacial width between the A and B species as shown in Fig. S6.We then determine χ by fitting these interfacial profiles using analytic expressions provided by Semenov (53): where ϕ A is the volume fraction of species A, z is the distance from the interface, and Here, R ee is the end-to-end distance of the neat homopolymer.This expression assumes two homopolymers of the same length N. We take the value of R ee = 29.1σLJ from neat homopolymer simulations with a bond angle rigidity κ = 2ϵ.This value is commensurate with the κ = 2ϵ simulations of Auhl et al. (55).We find excellent agreement between the analytical theory and our simulations.

THE EFFECT OF VARYING NON-BONDED CUTOFF RADIUS ON TOUGHNESS
In addition to the effect of entanglements on toughness, one can ask what effect our attractive component plays in transferring stress through entanglements and across load-bearing strands.To investigate, we take equilibrated, glassy configurations for a select set of simulations and modify their Lennard-Jones cutoff radius from R c = 2 .5σLJ to R c = 1.5σLJ .We then strain these simulations to failure using the same protocol as described previously and compare their toughness to that of the R c = 2.5σ LJ simulations in Fig. S7.These simulations do not suggest a systematic deviation from a perfectly correlated Γ R c =2.5σ LJ = Γ R c =1.5σ LJ relationship.These results suggest the attractive component in our simulations has little bearing on stress transfer, indicating topological constraints are much more important.

THE EFFECT OF THE BONDED POTENTIAL ON TOUGHNESS
The Morse potential is a convenient bonding potential that better approximates the potential energy of a diatomic molecule than harmonic bonds.As such, one can ask how changing from a FENE to a Morse-like potential would change the results of our work.Analogous to the previous section, we take glassy configurations for a select set of simulations and modify their bonding potential to a Morse potential before strain.Specifically, we used U b (r) = D 1 − e −α(r−r 0 ) 2 , where D = 17.3ϵ, α = 10.0/σLJ , and r 0 = 0.961σ LJ .This potential has the same minimum as the FENE potential and the same energy difference at r = 1.2σLJ , at which point the bond permanently breaks.We strain these simulations to failure and compare their toughness to that of the FENE bonded strain simulations in Fig. S8.Here, we find a small drop in the toughness of the compatibilized blends.Still, overall, the toughness of simulations with this bonded potential strongly correlates with the FENE-bonded simulations.These results suggest that switching to a Morse-like potential may quantitatively change parameters such as ϵ but will not modify our theory.
Fig. S1.SCFT parameterization assumptions.(A) This figure plots the change in the entanglement length scale as a function of the distance away from the interface (N e (z)) normalized by the bulk N e at Σ = 0.045σ L − J 2 and χN = 650.The block lengths L are listed in monomers.These results depend on Σ, χ, and compatibilizer architecture.Error bars represent standard errors of 3 replicates.(B) This image shows entanglement correlations calculated from Eq. S1 for compatibilizers with a block length of L = 32 monomers at Σ = 0.045σ L − J 2 and χN = 650.
Fig. S4.Validation of the use of load-bearing entanglements.This plot shows the true toughness against predicted toughness for a set of polymer blends compatibilized by short diblocks in which Σ is varied.Black circles show the predictions when using load-bearing entanglements in the model.Grey triangles show the same results when reanalyzed with all entanglements.The black line is the fit shown in Fig.3.
Fig. S6.Interfacial profiles.Plots of the interfacial profiles of a set of simulations without compatibilizers in which ϵ AB LJ = 0.5 (blue pentagons), 0.97 (pink squares), and 0.985ϵ LJ (red triangles).Points correspond to symbols used in Fig. 3 in the Main Text.The lines correspond to fits to analytic expressions by Semenov (53).

Fig. S7 .
Fig. S7.Toughness comparison for R c = 2.5σ LJ and R c = 1.5σLJ .Image plots toughnesses for simulations for block lengths L = 32, L = 64, and L = 128 monomers for our (N = 256, χN = 650) polymers.Before deformation, we modify the cutoff of the nonbonded Lennard-Jones potential from R c = 2.5σ LJ to R c = 1.5σLJ and compare their toughnesses, Γ R c =2.5σ LJ and Γ R c =1.5σ LJ .The line represents a perfect Γ R c =2.5σ LJ = Γ R c =1.5σ LJ relationship.We see no systematic deviation from this relationship.